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Abstract 

We investigate the convergence towards periodic orbits in discrete dynamical systems. We examine the probability that a 
randomly chosen point converges to a particular neighborhood of a periodic orbit in a fixed number of iterations, and we 
use linearized equations to examine the evolution near that neighborhood. The underlying idea is that points of stable 
periodic orbit are associated with intervals. We state and prove a theorem that details what regions of phase space are 
mapped into these intervals (once they are known) and how many iterations are required to get there. We also construct 
algorithms that allow our theoretical results to be implemented successfully in practice. 
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Introduction 

Periodic orbits are the most basic oscillations of nonlinear 
systems, and they also underlie extraordinarily complicated 
recurrent dynamics such as chaos [1-5]. Moreover, they occur 
ubiquitously in applications throughout the sciences and engi- 
neering. It is thus important to develop a deep understanding of 
periodic dynamics. 

It is important and common to question how long it takes a 
point in phase space to reach a stable periodic orbit from an 
arbitrary initial condition. When studying synchronization and 
other forms of collective behavior, it is crucial to examine not only 
the existence of stable periodic orbits but also the time that it takes 
to converge to such dynamics in both natural and human-designed 
systems [6-8]. For example, it is desirable to know how long it wUl 
take an engineered system that starts from an arbitrary initial 
condition to achieve the regular motion at which it is designed to 
work [9,10]. A system can also be perturbed from regular motion 
by accident, and it is important to estimate how long it will take to 
return to regular dynamics. Similar questions arise in physics 
[11,12], biology [6,13,14], and many other areas. It is also 
important to consider the time to synchronize networks [15-17] 
and to examine the convergence properties of algorithms for 
finding periodic orbits [2,18]. 

To study the problem of convergence time to periodic orbits, 
let's first consider the Hartman-Grobman Theorem [19,20], 
which states that the flow of a dynamical system (i.e., a vector 
field) near a hyperbolic equilibrium point is topologicaUy 
equivalent to the flow of its linearization near this equilibrium 
point. If all of the eigenvalues of the Jacobian matrix evaluated at 
an equilibrium have negative real parts, then this equilibrium 
point is reached exponentially fast when one is in a small 
neighborhood of it. To determine convergence time to a 
hyperbohc equilibrium, we thus need to calculate how long it 



takes to reach a neighborhood of the equilibrium from an arbitrary 
initial condition. After reaching the neighborhood, the temporal 
evolution is then governed by a linear dynamical system (which 
can be solved in closed form). An analogous result holds for 
hyperbohc periodic orbits in vector fields [21]. To turn periodic 
orbits in vector fields into fixed points in maps, one can use 
Poincare return maps, which faithfully capture properties of 
periodic orbits. A Poincare map can be interpreted as a discrete 
dynamical system, so the problem of determining how long it takes 
to reach a hyperbolic stable periodic orbit from arbitrary initial 
conditions in a vector field is reduced to the problem of 
determining how long it takes to reach the neighborhood of a 
hyperbolic fixed point in a discrete dynamical system. 

Our work considers how long it takes to reach a periodic orbit of 
a differential equation — starting from an arbitrary point in phase 
space — by using a Poincare return map of its associated vector 
field. For simplicity, suppose that a return map (which is built from 
a Poincare section) is unimodal. If we approximate the unimodal 
Poincare map by using a unimodal function f(x), then we can use 
f(x) in our algorithm to estimate the convergence time to the 
periodic orbit. Periodic motion is ubiquitous in models (and in 
nature), and it is important to explore how long it takes to 
converge to such behavior. 

In this paper, we prove a theorem for the rate of convergence to 
stable periodic orbits in discrete dynamical systems. Our basic 
strategy is as follows. We define the neighborhood Ip of a 
hyperbolic fixed point, and we calculate what fraction oj of the 
entire phase space / is mapped into Ip after q iterations. Using 
fi{w) and fi(I), respectively, to denote the measures of w and /, a 
point that is selected uniformly at random from / has a probability 
of ii(w)/ to reach Ip in q iterations. To illustrate our ideas, we 
will work with a one-dimensional (ID) discrete dynamical system 
x„ + 1 = f(Xn ; r) that is governed by a unimodal function / and is 
parametrized by a real number r. We focus on unimodal functions 
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for two primary reasons: (i) many important results in dynamical 
systems are based on such functions; and (ii) it is simpler to 
illustrate the salient ideas using them than with more complicated 
functions. 

To determine the set that is mapped into Ip, we take advantage 
of the fact that points in periodic orbits are repeated periodically, 
so their corresponding neighborhoods must also repeat periodi- 
cally. In theory, an alternative procedure would be to iterate 
backwards from Ip, but this does not work because one cannot 
control successive iterations of/^'. The function / is unimodal, so 
it is not bijective and in general one obtains multiple sets for each 
backward iteration of a single set. The number of sets grows 
geometrically, and one cannot in general locate them because an 
analytical expression for is not usually available. 

To explain the main ideas of this paper and for the sake of 
simplicity, consider a stable periodic orbit Op of period p that is 
born in p saddle-node bifurcations of /''. Every point x, (with 
(e{l, ...,/?}) of Op has a sibling point x* that is born in the same 
saddle-node bifurcation. Because f^{xj) = Xj and f^(x*) = x*, it 
follows that f''{Ii) = Ii, where /, = That is, x,, x], and /,■ all 

repeat periodically. Roughly speaking, we will build the interval Ip 
from the interval /,-. 

Consider a plot in which points along the horizontal axis are 
mapped via / to points along the vertical axis (as is usual for ID 
maps). The orbit Op is periodic with period p, so XjeOp implies 
that ^ = 0,1, . . .} yields p periodic points with a 

horizontal axis location oi Xj. We say that these points are located 
in the "column" Xj. Because f''{Xj) = XiSOp for some q, we obtain 
p points located in the same column Xj. These points are given by 
{(x//*(xy)) = (x/,x,), f = 1, . . . ,/)}. As we have indicated above, 
each point (xy,x,) is associated with an interval /,. No matter how 
many iterations we do, the fact that the orbit is periodic guarantees 
that there are exactly p intervals in the same place (where the 
points {xj,Xi) are located). We thereby know the exact number and 
locations of all inter\'als. 

To complete the picture, we must also take into account that if 
there exists an interval Wq^j such that f''{Wq,j) = Ii, then any point 
of Wq,i win reach a point of /,• in at most q iterations. The 
geometric construction above yields the inter\'al Wq^^, as one can 
see by drawing a pair of parallel line segments that intersect both 
f and the endpoLnts of the interval We will approximate /* by 
a set of such line segments so that we can easily calculate the 
intersection points. 

The remainder of this paper is organized as follows. First, we 
give definitions and their motivation. We then prove theorems that 
indicate how long it takes to reach the interval /, from an arbitrary 
initial condition. We then construct algorithms to implement the 
results of the theorems. Finally, we discuss a numerical example 
and then conclude. 

Definitions 

Consider the discrete dynamical system 

x„ + \=f{x„\r), f-.I^I, I=[a,b], (1) 

where f(x; r) is a one-parameter family of unimodal functions with 
negative Schwarzian derivative and a critical point at x=C. 
Without loss of generality, we suppose that there is a (both local 
and global) maximum at C. At a critical point of a map /, either 
/ =0 (as in the logistic map) or / does not exist (as in the tent 
map). Some of the results of this paper related with critical points 



only require continuous functions, which is a much weaker 
condition than the requirement of a negative Schwarzian 
derivative. 

Remark 1 Because f has a negative Schwarzian derivative, f does as 
well (because it is a composition of functions with negative Schwarzian 
derivatives). By using the chain rule, we obtain /''=0 only at externa. 
Therefore, # 0 between consecutive extrema. The Minimum Principle [22] 
for a function with negative Schwarzian derivative then guarantees that there is 
only one point of inflection between two consecutive extrema off^. If there were 
more than one point of infection, thenf^ =0 at least two points. One of them 
would be a maximum of j'^ , and the other one would he a minimum. This 
contradicts the Minimum Principle. Consequently, the graph of f between two 
consecutive externa has a sigmoidal shape (i.e., it looks like \ or\), which 
becomes increasingly steep as q becomes larger. This fact makes it possible to 
approximate f between two consecutive extrema by a line segment near the 
only point of irflection that is located between two consecutive extrema. 

Because the Schwarzian derivative of / is negative. Singer's 
Theorem [23] ensures that the system (1) has no more than one 
stable orbit for every fixed value of the parameter r. Additionally, 
the system (1) exhibits the well-known Feigenbaum cascade [24- 
26], which we show in the bifurcation diagram in Fig. 1. 

For a particular value of the parameter r, the map has p 
simultaneous saddle-node (SN) bifurcations, which result in an SN 
/7-periodic orbit. As r is varied, the SN orbit bifurcates into a stable 
orbit {SjYi^ 1 and an unstable orbit { UiYi^ l ■ "^^'^ points 5, and Ui 
are, respectively, the node and the saddle generated in an SN 
bifurcation, so Ut is the nearest unstable point to S,- (see Fig. 2). In 
other words, the points in the stable orbits (called "node orbits") 
are node points, whereas the points in the unstable orbit (called 
"saddle orbits") are saddle points. From Remark 1, we know that 
the neighborhoods of these points are concave or convex. 

The derivative of f is 1 at the fixed point where the SN 
bifurcation takes place. As one varies r, the derivative evaluated at 
that bifurcation point changes continuously from I to — 1 . When 
the derivative is —1, the stable orbit (i.e., the node orbit) 
undergoes a period-doubling bifurcation. As a result, the stable 
orbit becomes unstable (yielding the orbit {J7,}f^[) and two new 
stable orbits ({•S,-, and {■S'/jIfj^i) appear. The points and 
5,2 are nodes, and the point f7, is a saddle. From our geometric 
approach, the intervals (5,-, , [/,) and {Si-, , [/,) that are generated via 
the period-doubling bifurcation behave in the same way as the 
interval (S,-,{7,-) that was generated in the SN bifurcation. 
Therefore, we can drop the indices "1" and "2" and write 
{Si, Ui) for the orbits that arise from both the SN bifurcation and 
the period-doubling bifurcation. 

Notation 1 Let Uf denote the nearest point to Si that results from the 
intersection of the line x„+i = Ui with f (see Figs. 2, 3, and 4). 

Definition 1 Consider the points x,- and x- that satisfji f^{Xi) = Ui 
and f^{Xi)= Uj. Iff is concave (respectively, convex) in a neighborhood of 
Si, we say that Ip^ = (x/,x,) [respectively, Ip. = (Xi,Xi)] is the ith capture 
interval of the stable p-periodic orbit {S,}^^ j and that Ip=\^ Ipt is the 

i 

(aggregate) capture interval of the stable p-periodic orbit {Si}^^^. 
Notation 2 Let Ip,^ denote the subinterval Ip, that contains the critical 

point C. 

From Definition 1, we see for all xe/p, that f"''{x)elpi and 
f"^(x)^Si as n^co. Iterations of points xe/p, are repelled from 
Ui and Uf, and they are attracted to S,-. The system (1) is 
linearizable around the fixed points Si and Ui. (Observe that 
f(Ui) = Ui, so we also have control over this point.) Consequently, 
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Figure 1. Bifurcation diagram of a unimodal map with a negative ScKiwarzian derivative. There is a period-doubling cascade on the left, 
and there are also period-doubling cascades inside several windows (the broad, clear bands) of periodic behavior. Saddle-node orbits arise at the 
onset of such windows in the chaotic area. 
doi:1 0.1 371/journal.pone.0092652.g001 



the convergence of iterations of xelp, to 5,- is governed by the 
eigenvalues of the Jacobian matrix Df. 

Because we can control the evolution inside /p, , we can examine 
how long it takes to reach Ip^ starting from an arbitrary point xel. 
As we will see below, to obtain this result, we need to discern 
which subintervals of / are mapped by into Ip. for arbitrary q. 
The first step in this goal is to split the interval / in which f is 
defined into subintervals in which is monotonic. 

Definition 1LetA = {q2,q3, ■ ■ . ,qk — i\q2< ■ • ■ <1k — ^} be the 
set of points at which f has extrema. Let B = {qi = a,qk = h}, and we 



recall that we are considering the interval I = [a,b] . We will call 

Pmon-r=A[JB={quq2,q3,---,qk-U , 

the partition ot monoto- 

qk\q\<qi< ■ ■ ■ <qk-\<qk} 

nicity of f . We will call I^j = [qj,qj + 1] (where j =\, . . . ,k—\) the jlh 
interval of monoto nicity of /*. 

By construction, / = {a,b\ = [_J/ij , and is monotonic in /^^,. As 

j 

we will explain below, one can calculate intervals of monotonicity 
Iqj easily by using Lemmas 1 and 2. 
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Figure 2. Geometric calculation of the three subintervals corresponding to a 3-periodic orbit (in blue). See Fig. 3 for a better view 
of the orbit. These subintervals are determined by the three pairs of black, horizontal, parallel line segments that intersect /'', U,, and f. (We only 
indicate one f/, in the figure.) One needs to take into account the intersection points of /'' with all 6 parallel line segments. See Figs. 4 and 5 for more 
detail. The plot in this figure uses the logistic map. The blue orbit is a period-3 supercycle and f*3. 83187405528331556841. 
doi:1 0.1 371/journal.pone.0092652.g002 
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Figure 3. The 3-periodic orbit f =f^ from Fig. 2. 

doi:10.1371/journal.pone.0092652.g003 
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Figure 4. Magnification of Fig. 2. We show the interval W^., in which/'' does not have any extrema in the region between the horizontal parallel 
lines. The horizontal line that crosses Uj and intersects with /'' determines (7,.. The vertical lines that intersect (7, and C/, determine .y, and x^, 
respectively. We obtain locations for the points .y,,,,^, and x^^^, because their images under the map/'' are x, and x,, respectively. We thereby construct 
the subinterval W,, . We depict the mapping of the subinterval W,i^^ using a filled green arrow the mapping of another subinterval using the filled 
blue arrow. 

doi:10.1371/journal.pone.0092652.g004 
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Figure 5. Another magnification of Fig. 2. We show the interval W^,,,, in which/'' has an extremum in the region between the horizontal parallel 
lines. The horizontal line that crosses and intersects /'' determines t/, . The vertical lines that intersect (7, and (/, determine the points x, and .y,, 
respectively. We obtain the locations for the points .t,,,,^, and x-,^. because their images under the map /'' are .y,- and x-, respectively. We thereby 
construct the subinterval Wq^^. We depict the mapping of the subinterval Wi^^^ using a filled green arrow. 
doi:1 0.1 371 /journal.pone.0092652.g005 



Once we know the intervals in which is monotonic, it is easy 
to obtain subintervals of / that are mapped by /' into /p,. 
We proceed geometrically (see Figs. 2, 4, and 5): 

• draw parallel lines through the points x,- and -X, (i.e., through 
the endpoints of /p,); 

• obtain the points at which the lines intersect 

• calculate which points are mapped by into the intersection 
points of (ii), for which one uses the fact that /' is monotonic in 

• determine, using the points obtained in (iii), the interval that is 
mapped by /' into //•,. 

Using this geometric perspective, we make the following 
definitions. 

Definition 3 Let 



(ii) Iff has extrema in Aij (see Figs. 2 and 5), then we let 



iXi,qj,Xi,qj), if X,,q.>X. 



(2) 



where the points Xi^q- ,Xi ^elq^ = [qj,qjj^\\, and they satisfy f''(x^qj) = Xi 
andfl{x[^,) = x]. 

(i) Iff does not have extrema in Aij (see Figs. 2 and 4), then we let 



(3) 



lij- 



iXi,qj,C\, 



'/ Xi^q. < C 

if x^. > C ' 



(4) 



Remark: If we did not take point (ii) into account, then f 
would not be monotonic in W^,, . 

By construction, all points xeWq^j reach /p, in at most q 

iterations (see Figs. 2, 4, and 5). That is, /'(fy^.,) = /p, for l<q. 

Definition 4 We call Wr, = [jWq,j the (/-capture interval of /p, , 
,/" 

as Ip^ is captured after at most q iterations. The interval 
Wr = Wji^ = Wq^^ is then the (/-capture interval of the orbit 

Observe that Wq^j can be the empty set for some values of j. 

Theorems 

Once we know Wr, we can calculate the probability that a 
point picked uniformly at random from phase space is located in 
Wr. We can then calculate the probability that that point reaches 
a capture interval of Op in at most q iterations. We let ii(Wr) 
denote the measure of Wr, and we have the following theorem. 

Theorem 1 Let 0^ = {5',}^^i be a stable p-periodic orbit of the 
.system (1). Given an arbitrary point xel, the probability to reach a capture 
interval of Op after at most q iterations is 
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nil) b — a 



(5) 



Proof 1 From the definition (4) oJWr, all xgWr satisfy f'{x)eIpfor 
l<q. There always exist values qfl<q such that f\x)elp because extrema 
of f'{x) that satisfy l<q are necessarily also extrema of and points 
belonging to the latter set of extrema reach a capture interval of Op after at mast 

I iteratiims (see Lemma 1 heloic). Consequently, one reaches I p from xgWr 
after at most q iterations (cmd we note that it need not be exactly q iterations). 
Thus, the probability to reach Ip from an arbitrary point xel after at most q 
iterations (i.e., the probability that xe Wr) is 



li(WR) 



KWr) 

b — a 



(6) 



Corollary 1 With the hypotheses of Theorem 1, the probability to reach 
a capture interval of Op in exactly q iterations is 



This answers the question of how long it takes to reach a 
capture interval of a /;-periodic orbit from an arbitrary point. 
However, we also need to calculate fi( Wr). To do this, we need to 
understand the structure of Wq^j. As the following lemma indicates, 
some of these subintervals are located where /* is monotonic and 
others contain extrema of f . 

Lemma I Iff : I->I is an unimodal function with a critical point 
at C, then /'(x) has extrema 

(i) at points for which /*^'(x) = C; 

(ii) at tlie same points at which '(x) has extrema. 
Proof 2 

(i) For all xel .such that f''^^(x)=C, we know that 
/*(x) ^'(x))=/(C). Therefore, has an extremum because f 
has an extremum. 

(ii) Write / = //, U {C} KJJr, where Jl = (a,C) and Jr = (C,b), so f 
is a monotonic fimctim on the intervals Jl and Jr. 

(ii.a) If xeJi or xgJr and the function f^^ix) has an extremum, then 
we know that f ' (x) is a monotonically increa.nng function on one .nde of X 
and a monotonically decreasing function on the other, (jmsequently , 
/*(x) ^ '(x)) is the compo.sition of two monotonic functions (f and 

f^^), both of which are increasing (or decreasing) on one side of x. On the 
other .side of X, one of them is increasing and the other is decreasing. Therefore, 
there is an extremum at x. 

(ii.b) Otherwise, if /'~'(x) = C, then we see straightforwardly that /' 
has an extremum. 

We have just seen how to determine the locations of extrema of 
/'.We also need to know the values that f takes at these extrema. 

As we will see below, if the system (1) has a stable /i-periodic 
orbit and q>p, then the values that f takes at its extrema are the 
same as those that f takes at its extrema. This makes it possible to 
calculate the subintervals Wgj that are associated with extrema of 
/' by using Ip.^ and the derivative of /. 

Lemma 2 Let Op be a stable p-periodic orbit of the system (1). The 
coordinates if the extrema of f (where q>p) are (x,c/'~'|^(C))_, where 
Xic denotes the points xe/ .such that f'(x) = C, the index i takes values of 
i = Q,\, . . . ,q—\ (where we note that f^ = Id is the identity map), and 
q — i\p = {q — i)mod p. 



Proof 3 According to Lemma 1, the extrema of f are 

(i) xel such that f''^\x) = C; 

(ii) xel such that /*^'(x) is an extremum. 
It thus follows that the extrema off^~^ are 

(ia) xel .such that f1^^(x)=C; 

(iib) xel such that f ^(x) is an extremum. 

Repeating the process, we obtain that extrema of f are located at Xic, 
where ! = 0,1, 1. The value of f at Xic is 

fHx)=p-'(f'{x))=f''-(C). 

Because Op is a .stahh' p-periodic orbit, there exists one point of Op near C 
that is repated periodically after p iterations. Consequently, 
{Cf{C)f^(C),...f~^(C)} is a periodic sequence and 
/"-'{€)= (O. 

Algorithms 

As we discussed above, Lemmas 1 and 2 determine intervals of 
mono tonicity (see definition 2), and they also make it possible to 
construct algorithms for calculating Wq^^. 

For these algorithms, we approximate f by line segments in the 
subintervals in which f is monotonic. This approximation is very 
good unless one is extremely close to an extremum (see Fig. 6), and 
this is already the case even for relatively small q (as we will 
demonstrate below). Additionally, recall that Wq^j is determined by 
the intersection points of f with line segments. Therefore, once 
we have approximated f by a set of line segments, it is 
straightforward to calculate those intersection points. 

Algorithm 1 (Calculating coordinates for extrema of /*) 
Suppose that we know the coordinates of the extrema of f'^~^ . According to 
Lemma f the extrema of f are located at the points 

(i) xel such that /*^'(x) = C and /*^'(x) is not an extiemum; 

(ii) xel such that /*^'(x) is an extremum. 

We know the extrema in (ii) by hypothesis. To find the extrema in (i), we 
need to calculate the points xel that satisfy ' (x) = C. Because we know 
the coordinates of extema off^~^, we construct the lines that connect two 
consecutive extrema (see Fig. 7). Let x„^\= ax„ + h be the equation for .such 
a line. We solve ax„ -\-b = C to obtain a .seed that we can u.se in any of the 
many numerous numerical methods for obtaining roots of rwnlinear algebraic 
equations. Observe that f is monotonic in the interval in which the line 
x„-i-i =ax„+b is defined. This circumvents any problem that there might 
otherwise be in obtaining a good seed to ensure convergence of the root solver. 
Moreover, we have as many seeds as there are points xel that satisfy 
./"'^ ' (x) = C. Mote that we need to construct both the line that connects 
{af\a)) with the first extremum off^^ and the line that connects {bf{b)) 
with the last extremum qff''~^. 

To calculate the points .xel for which /"'^'(x) is an extremum, we apply 
this algorithm recursively, and we note that we know by hypothesis that f has 
an extrememum at C. We first build the line segments that connect {af(a)) 
with {Cf{C)) and {CfiC}) with {bf{b)). These two line segments give 
seeds from which to determine the points xel that satisfy f(x) = C. We 
thereby obtain the coordinates for the extrema of f^. We then use the same 
procedure to obtain the coordinates for extrema of f^ ,f*, ■ ■ ■ ,/*. 

We will see below that if the system (1) has a stable ^-periodic 
orbit and q»p, then the points xel with /^^'(x) = C are given to 
a very good approximation by the intersection points of two lines. 
Moreover, as one can see in Fig. 6, the value ^ = 6 is already large 
enough to approximate f very successfully by a set of line 
segments when / is the logistic map. 

Algorithm 2 (Calculation of Wq.^ in the system (1)) Suppose that 
we know the coordinates of the extrema of /* (e.g., by computing them using 
Algorithm 1). We want to obtain Wq-from the definition (3), where 
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Figure 6. Outside of tKie intervals W^.., we approximate tKie map f using line segments. A line segments connects the upper endpoint of 
the interval W^^. to the lower endpoint of Wq^.^^. The map is very well approximated using line segments as long as one is not too close to an 
extremum. We again use the logistic map to illustrate our procedure. The blue curve is a period-6 supercycle and ;-3s3. 99758311825456726610. See 
Fig. 7 for a magnification of this figure. 
doi:10.1371/journal.pone.0092652.g006 



/*fe,,P = x,-, /''(x,„,) = .v-, 



.9,' 



(7) 



Taking into account that Xic is known, we construct the sequence 
Sic = {xifl,Xi i , . . . ,Xij = C} , 



and the ith capture interval is 

(Xi,x'i), if 
(x[,Xi), if 



Xi<X^, 
Xi>X^. 



(8) 



To determine the points Xj ^j and x-^ , we first approximated them by 
replacing f by line segments that connect consecutive extrema of f (i.e., by the 
same procedure that we use in Algorithm 1 to obtain approximations of points). 
Using the approximations of Xi ^j and x^^^, we construct the interval 
lapp = (Xi.qj yX^ and then check if there is an extremum of in /app- (This 
is trivial because we know the coordinates of the extrema of f .) We need to 
consider two cases. 

(i) The map has no extrema in /app- Tliis is equivalent to case (i) of 
Algorithm 1. We use the approximations of Xi q^ and x^^, as seeds in a 
numerical root-finding method. 

(ii) The map f has extrema in /app- This is equivalent to ca.se (ii) of 
Algorithm 1. 

If there is an extremum of f in /app, then that extremum is necessarily one 
of the extrema given by Lemma 2: (x,c/* 'I^(C)). Because f'(xic) = C 
and f is a continuous function, there must exist an interval lic such that 
XicsIiC andf'(Iic)<^Ipic 



where Xic = Xifi, Xi^k=f''ixic), andf°(Xic) = Xic 
Let Li ic be the linear map whose graph is the line of slope f (x,;^) that 
intersects the point Xi^^. If the period p of the orbit is sufficiently large, then we 
can approximate f near xi^k (where fc = 0,l, — Ij by the linear map 
Li k. Thus, instead of iterating //c with the map f to obtain Ip,^., we iterate 
liC with the linear map Li^^ that approximates f. That is. 



Ipj (J « Lij- 1 . . . Lifi(Iic) ■ 



Because each Lijc is a linear map, it is .straightforward to compute Z,,- ^, 
and hence to compute 



liC^L-g . ..L. -liilp.^). 



At the end of this section, we will discuss the error that is introduced by this 
approximation. 
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Figure 7. Graphs of (blue) and (red) for the same value of the parameter r (when / is the logistic map) as in Fig. 6. The darl< pink 
line joins two consecutive extrema of and the black line is the tangent line that crosses through the inflection point. Both lines are approximations 
to/''. As expected, the approximation is better for the larger value of q. 
doi:10.1371/journal.pone.0092652.g007 



The interval lie tliat we have just constructed is the interval 



W, 



'lit ' 



{x^.,C\, if 
[C,x^.), if 



if < C 
if > C 



(9) 



Fig. 7). The Taylor polynomial of degree 3 of f around the 
inflection point Xm is 



fix) ^/«(Xi„f ) +/"' {X,,c){x - X,nf ) + ^Z"'" {X,„c){x - X,,c? ■ 



that we seek. 

In Algorithm 1 , we constructed line segments that connect two 
consecutive extrema of They are located in the intervals 
[qj,qj+\) and [qj+\>qj+i), respectively. We now have intervals 
Wii.j<^[qj,qj+\) and Wqij + \ <^[qj+\,qj+2) that contain these two 
consecutive extrema of f, so we construct the line segment that 
connects the upper endpoint of PF^,^. to the lower endpoint of 
Wjij+i . (Note that we do not connect the two extrema directly via a 
line segment.) For ^ » 1 , this line segment approximates f outside 
of the intervals Wq.j and W^^j^^. See Fig. 6, which illustrates (for 

the case when / is the logistic map) that we can approximate f^ by 
a set of line segments for q = 6. We can then use these line 
segments in Algorithm 1, and we do not need numerical 
computations to find the intersection points. 

As we discussed previously, we can replace f by linear 
expressions to approximate the intersection points when deter- 
mining Wr in Algorithms 1 and 2. Replacing /' by a linear 
approximation simplifies operations and reduces the amount of 
calculation. To determine the desired intersection points, we have 
thereby replaced a numerical method for obtaining roots of 
nonlinear algebraic equations by an analytical calculation that uses 
a system of two linear equations. We now estimate the error of 
replacing/* by lines segments. The line segments that replace the 
function f intersect /' very close to the unique point of inflection 
between a pair of consecutive extrema of /'' (see Remark 1 and 



Consequentiy, the error of approximating /* by the line 

/*(^inf ) +f'''(Xmt)(x - Xi„t) is 

Error = | ^/""'(xirf)(x-Xi„f)'| « | ^/*"'(xi„f) (^^^l, (10) 

where we have taken into account that there are more tiian 2' 
local extrema of f in the interval [a,h] . The exponential growth of 
2* enforces a fast decay in the error. Consequently, using line 
segments to approximate /' is an effective procedure with only a 
small error. 

Numerical Example 

Algorithms 1 and 2 are based on the same procedure: 
approximate fix) by a line y{x) = ax + h and solve y(x) = C to 
obtain an approximation of the f^ix) = C (instead of solving 
f''(x) = C directly). In this section, we consider an example 
application of Algorithm 1. 

To obtain the critical points of /* + 1 , we need to calculate the 
points that satisfy f = C. Suppose that q = 6 (and again see Fig. 6 
for an illustration of the line-segment approximation with q = () for 
the logistic map). The biggest distance between consecutive 
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extrema occurs near the critical point C, so we approximate / by 
a line segment in this region to obtain an upper bound for the 
error. The extrema are located at (4.525 x 10^',2.414 x 10^^) 
and (4.787 x 10~\9.994 x 10~'), and they are connected by the 
line v« 38.053 x—l6, from which we obtain the approximation 
•Xapp* 0.453 for the solution of/^(x) = C. From direct computa- 
tion, the value of x that satisfies f^(x) = C is 0.465. The 
relative error is £'rel~2.58%, and this is the largest error in this 
example from all of the approximating lines segments. As we 
showed in equation (10), the error decreases exponentially. Hence, 
when we approximate using line segments, the relative error 

win be bounded above by iirel *2.58/2™ %. One can observe this 
decrease in error in Fig. 7, in which we plot both and for the 
logistic map and the same parameter value r. Observe that several 
extrema of lie between onsecutive extrema of so using the 
hne-segment approximation in induces a much smaller error 
than using it in f^. 

Conclusions and Discussion 

When studying dynamical systems, it is important to consider 
not only whether one converges to periodic orbits but also how 
long it takes to do so. We show how to do this explicitly in one- 
dimensional discrete dynamical systems governed by unimodal 
functions. We obtain theoretical results on this convergence and 
develop practical algorithms to exploit them. These algorithms are 
both fast and simple, as they are linear procedures. One can also 
apply our results to multimodal one-dimensional maps by 
separately examining regions of parameter space near each local 
extremum. 



Although we have focused on periodic dynamics, the ideas that 
we have illustrated in this paper can also be helpful for trying to 
understand the dynamics of chaotic systems. Two important 

properties of a chaotic attractor are that (i) its skeleton can be 
constructed (via a "cycle expansion") by considering a set of 
infinitely many unstable periodic orbits; and (ii) small neighbor- 
hoods of the unstable orbits that constitute the skeleton are visited 
ergodically l)y dynamics that traverse the attractor [4]. In Refs. 
[21,22], Schmelcher and Diakonos developed a method to detect 
unstable periodic orbits of chaotic dynamical systems. They 
transformed the unstable periodic orbits into stable ones by using a 
universal set of linear transformations. One could use the results of 
the present paper after applying such transformations. Moreover, 
the smallest-period unstable periodic orbits tend to be the most 
important orbits for an attractor's skeleton [4], so our results 
should provide a practical tool that can be used to help gain 
insights on chaotic dynamics. 

Once unstable orbits has been transformed into stable ones we 
can use results of this paper to answer the above question. 
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